Geometric visualization of self-propulsion in a complex medium 
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Combining geometric mechanics theory, laboratory robotic experiment and numerical simulation, 
we study the locomotion in granular media (GM) of the simplest non-inertial swimmer, the Purcell 
three-link swimmer. Using granular resistive force laws as inputs, the theory relates translation and 
rotation of the body to shape changes (movements of the links). This allows analysis, visualization, 
and prediction of effective movements that are verified by experiment. The geometric approach also 
facilitates comparison between swimming in GM and in viscous fluids. 



Introduction-Locomotion of animals and robots 
emerges through the interplay of body deformations cou- 
pled to an environment. Finding this relationship is often 
a challenge: for example, in Newtonian fluids, although 
researchers have long analyzed [lHl] the Navier-Stokes 
equations and simpler representations Q to gain insight 
into flight and swimming, analytic investigation is often 
impossible, and high fidelity approximations are compu- 
tationally costly. Studying the motion of organisms [B[ 
and robots [6[ that maneuver through complex environ- 
ments like sand, rubble, and debris, and microscopic or- 
ganisms Q that move through complex biomaterials can 
be even more complicated — often such materials are not 
even described by equations at the level of Navier-Stokes. 

Certain kinds of movement are kinematic, in that the 
net displacement is a function of the deformation and is 
independent of its rate. Taking advantage of this prop- 
erty for low Reynolds number (Re) swimming in viscous 
Newtonian fluids, Shapere and Wilczek p| introduced a 
geometric approach using the notion of gauge symme- 
tries, which are equivalencies in the system dynamics 
across different configurations. These symmetries reduce 
the effective dimensionality of the system and facilitate 
interpretation of dynamics in terms of geometric concepts 
such as areas, lengths, and curvatures. The geometric 
approach has been further developed d, Q and now en- 
ables evaluations [T3| of systems' locomotion capabilities 
in the form of low-dimensional, readily visualizable rep- 
resentations of the system's motion for any gait (cyclic 
change in body shape). These tools allow useful gaits 
to be identified by inspection, without costly trial-and- 
error optimization However, the insights afforded by 
these geometric tools have been restricted to systems — 
including viscous swimmers 0, & El, El and planar reori- 
enting satellites 14j — with analytically-describable linear 
dynamics. 

Previously, we studied [HI the sand-swimming of a 
sandfish lizard and a robot model in dry granular media 
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FIG. 1. Three-link swimmers: (A) The robot resting on a 
granular medium (GM) of plastic spheres. The lit masts at 
the ends of the links are markers for tracking the robot when 
its body is submerged. The rubber skin on the front section 
has been removed to display the mechanical structure. Each 
link is 5.4 x 2.8 x 14.7 cm 3 . (B) Analytical three-link model 
with a body frame corresponding to a weighted average of the 
link positions and orientations, (c) Force relations for a rod 
dragged in GM (adapted from [15|]) [solid] and a long slender 
rod dragged in a low Re fluid [dashed]. Scaling is chosen such 
that the maximum values of Fu are equal for the GM and low 
Re systems. 



(GM), arguably the simplest flowing terrestrial material. 
In theory, simulation and experiment, we demonstrated 
that the GM surrounding an immersed undulatory swim- 
mer could be modeled as a "frictional fluid" in which 
forces are dominated by Coulomb friction, making them 
insensitive to rate and in which inertial effects are small. 
However, our analysis could only hint at the range of be- 
haviors possible in sand-swimming. Here we demonstrate 
the efficacy of the geometric approach to reveal princi- 
ples of swimming in GM, despite a lack of fundamental 
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equations of motion. We empirically generate a geomet- 
ric swimming model in GM for the three-link swimmer 
(Fig. [TJ hrst introduced by Purcell as a simple swim- 
mer to study locomotion in viscous fluids [l7|- We use 
this model to analyze different locomotor behaviors and 
to compare swimming in GM to swimming in a viscous 
fluid. 

Geometric mechanics, the resistive force model, and 
the three link swimmer-The key ingredient in applying 
geometric theory to motion in GM is the ansatz that 
at any given shape (with joint angles specified by the 
vector a — (a.\,ot2)), the swimmer's body velocity £ is 
linearly proportional to its shape velocity (at), such that 
the relationship between shape, shape velocity, and body 
velocity can be expressed as 

£ = A(a) • at = (A x (a) ■ a, A v (a) ■ a, A e (a) ■ at), (I) 

where A(a) is referred to as the local connection (or Ja- 
cobian) matrix 18]. Local connection models have been 




identified for diverse locomotion modes 
l2lT ] , in particular for swimmers in low Re 7 , 
viscous environments that qualitatively resemble those 
seen in granular swimming |15| . The existence of such a 
model for motion in GM is further suggested by our pre- 
vious results [HI showing that sand-swimming is kine- 
matic. Local linearity between shape and position veloc- 
ities sufficies to produce kinematic motion, and is only 
slightly stronger than the necessary condition, propor- 
tionality between body and shape velocities. 

GM lack equations equivalent to Navier- Stokes, so an- 
alytic derivations of the local connection used previously 
are not applicable. In their place, we have developed a 
numerical means of identifying A, based on our empir- 
ically obtained granular resistive force laws [1, [H| and 
the observation that inertial forces on a low-speed swim- 
mer are sufficiently small that the swimmer moves quasi- 
statically. The force laws (Fig. [IJ resemble those in low 
Re fluids [22} , although forces perpendicular to body seg- 
ments are enhanced relative to those in true fluids. In- 
tegrating these forces along the swimmer's body at dif- 
ferent a, a, and £ combinations and solving for the com- 
binations that give force equilibria yields mappings from 
shape velocity to body velocity at each shape. 

We find that across the space of shapes, graphs of the 
components of £ as functions of a are nearly planar, as 
in Fig. [2j This planarity means that, despite the lack of 
an analytic linear relationship between shape and body 
velocities, there exists an implicit linear relationship to 
which we can fit a local connection of the form in Eqn. 1. 
Each row of A (a) is then obtained by finding the best- fit 
planes for each component of the body velocity at each 
shape a [23[. Taking these components over the set of 
shapes, we visualize each component of A as a vector 
field on the shape space (Fig. 3a), termed a connection 
vector field. 
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FIG. 2. The components of £ for a particular initial shape 
(inset), with a — (0.85, —1.14) radians; velocities are given 
in terms of body-lengths (£ x and £ y ) or radians (£ t heta) per 
second. The components of the equilibrium £ (shown here for 
a body frame attached to the middle link) are almost linear 
functions of the shape velocity, and so can be characterized 
by the gradient vectors of the best-fit planes (black arrows). 
The color map shows the agreement between the equilibrium 
£ and the fitted plane, with light regions exhibiting the least 
error. 
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FIG. 3. (A) The local connection vector fields (Eqn. 1) 
A. x ,A. y ,A. . The red path through the shape space denotes 
a circle gait, and the grey circle identifies the gradients taken 
from Fig. 2. (B) Sequence of shapes and displacements along 
the gait shown above. The vertical red line is a reference po- 
sition for displacement. (C) Along the gait, the £ and net dis- 
placement predicted by the linear model (thick red line), the 
DEM simulation (dotted line) and robot experiments (thin 
black line) for a representative experiment at frequency 0.5 
Hz. Lengths and angles are given in terms of body-lengths 
and radians, respectively. Time is in fractions of a gait cycle. 



Testing the linear theory- To test if the linear approx- 
imation (local connection approach) of the system dy- 
namics can accurately compute movement over different 
gaits, we constructed a three-link robot in both experi- 
ment and multi-particle discrete element method (DEM) 
simulations, using identical techniques and parameters to 
those reported in [6[. The robot (Fig. 1) consists of three 
wooden segments connected by two servomotors (Hitec, 
HSR 5980SG); the total mass of the robot is 0.56 kg. The 
segments are covered by a thin latex sleeve giving the 
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robot a body-particle coefficient of friction of 0.4. The 
robot is fully immersed in a large bed of 5.87 ± 0.06 mm 
diameter plastic spheres at a depth of 5.5 cm from the 
top of the robot. Motion was tracked through a camera 
via the position of LED masts. Data was collected and 
averaged over 3 runs. The DEM simulation used 6 mm di- 
ameter particles and a particle-particle and particle-body 
collision model incorporating Hertzian contact, normal 
dissipation and tangential Coulomb friction, with param- 
eters previously validated against experiment 

Si- 
Integrating the linear model's prediction of the swim- 
mer's motion during a shape change is equivalent to tak- 
ing line integrals on the connection vector fields along 
paths the system traces. For example, consider a gait 
in which the joints oscillate sinusoidally with a quarter- 
phase offset, producing a traveling wave of deformation 
along the body. This gait traces out a circle in the space 
of joint angles, shown in Fig. [3H, generating positive 
and negative £ as it flows along and against the vector 
fields, depicted in Fig. [3]C; these velocities can then be 
integrated into net displacements relative to the start- 
ing body frame. As illustrated in Fig. [3]C, the x and 
9 components of the velocities and integrated positions 
predicted by the linear model agree with those found in 
experiment and DEM. The y component of the veloc- 
ity and displacement from the linear model qualitatively 
agree with those from the experiment and DEM but the 
magnitudes are reduced. We will return to this point. 

Constraint curvature functions-The local connection 
model of granular swimming simplifies evaluation of the 
displacement produced by a gait, but its greater useful- 
ness derives from the ability to calculate the associated 
constraint curvature functions (CCFs). These functions 
(calculated following the procedure in 0, [2^ ) are plot- 
ted for GM in Fig. gJ4. CCFs are closely related to the 
curls of the connection vector fields. In an extension 
of Stokes theorem, this curl-like nature means that the 
net displacement induced by a gait corresponds to the 
area integrals of the CCFs over the region the gait en- 
closes in the shape space. Stokes theorem only applies 
for systems where the integrations are commutative and 
so does not directly apply to our swimmers, for which 
body-frame translations and rotations do not commute. 
CCFs circumvent this limitation by augmenting the curl 
with a Lie bracket term that linearly approximates the 
effects of noncommutativity. We introduced [l(| coordi- 
nate optimization techniques that minimize the error in 
this linearization; for the granular-swimmer considered 
here (as in [l(| and [13]) the error is negligible. 

By visualizing gaits as enclosures of area, the CCFs 
provide a comprehensive overview of how gait patterns 
interact with system constraints to produce net displace- 
ment. For example, we can explain the net forward dis- 
placement for the gait in Fig. [3]a by overlaying it on the 
x b CCF from Fig. gH (see inset of Fig. The large 
negative region at the center of the x b plot indicates that 




FIG. 4. Fundamental geometric diagrams, the constraint cur- 
vature functions (CCFs) for the (A) granular three- link swim- 
mer and (B) low Re swimmer of equivalent dimensions (bot- 
tom). Each swimmer's body frame is optimized to a weighted- 
average of the three link frames, see [2^]. The units of the 
field values are body-lengths/joint-radians 2 for x b and y b and 
orientation-radians/joint-radians 2 for 9 b '. Values on the x b 
and y b fields have been multiplied by 100. 



cycles in this region produce net x translation relative to 
the starting frame, with clockwise (negatively oriented) 
cycles generating positive displacement. See SI Movie 1. 
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FIG. 5. Predictions of movement for the three-link swimmer: 
(A) CCF estimate of displacement compared to experimental 
and DEM results as a function of stroke amplitude (radius 
of the circle traced in the shape space). Dashed horizontal 
lines indicate displacements for the butterfly gait in DEM and 
RFT. The inset shows the circle gait overlaid on the x b CCF, 
along with a graphically-optimized "butterfly" gaits (dashed 
path). (B) CCF estimate of net rotation for figure-eight gait 
of different stroke amplitudes, the radius of one circle of the 
figure-eight overlaid on the 8 b CCF in the inset. Frequency 
was 0.5 Hz in (A) and 0.17 Hz sin (B). 

The CCFs also give insight into the relationship be- 
tween net displacement and magnitude of the joint mo- 
tions. For small amplitude circular gaits, where the sign 
is negative, the net displacement scales approximately 
quadratically with amplitude, tracking the rise in the en- 
closed area. At large amplitudes, the gait includes pos- 
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itive regions near the corners of the plot, reducing the 
area integral. These behaviors point to the optimal am- 
plitude as that which encompasses as much of the central 
negative region as possible while avoiding the outlying 
positive regions. This geometric interpretation also sug- 
gests gaits, such as the butterfly in the inset of Fig. [5]A, 
that produce more displacement than any circle by better 
conforming to the sign-definite regions; such a gait was 
created by fitting a polynomial curve to a set of points 
that were on the zero-set, symmetric about the origin, 
and avoided extreme joint angles [26[ . 

Agreement between theory, experiment and simulation 
is excellent at small amplitudes and good at larger am- 
plitudes. However, while the CCFs qualitatively predict 
forward movement at the highest amplitudes, quantita- 
tive agreement is lacking. Previous work [l(J Ell suggests 
a cause for this error. The velocity error is largest when 
joint angles are near the ai — cii line and away from 
the origin, i.e., with links in a "C" shape. This error 
appears in Fig. [3]C as the system experiencing greater- 
than-expected velocity in the body x and y directions in 
the vicinity of this line (near t = 0.5). In similar sys- 
tems (e.g. the "kinematic snake" 1(J III), this line is a 
kind of kinematic singularity, which forces constraints to 
slip. For granular-swimming, in which forces are inher- 
ently nonlinear, we postulate that the singularity makes 
the system reach different equilibria. 

The CCFs facilitate study of movements that are rela- 
tively unexplored in swimming locomotion, for example, 
turning in place, see SI Movie 2. To design this gait, we 
observe that the x component of the CCF is evenly sym- 
metric around both the ol\ = ai and oi\ = —cti lines, the 
y component is odd around u\ — ai and even around 
u\ = — a>2, and the 9 component is even around a.\ = «2, 
but odd around ct\ — —ot2. The figure-eight gait depicted 
in the inset of Fig. [5]S will therefore produce a net ro- 
tation of the system — the two loops are in opposite di- 
rections in oppositely-signed regions, and so their effects 
add in 9 while canceling in x and y. In this gait, agree- 
ment between theory and experiment/DEM is excellent; 
we hypothesize that here the singularity described above 
is avoided. 

The CCFs also facilitate comparison of swimming in 
different environments. The granular swimmer's CCFs in 
Fig. EH share a structure with those of other three-link 
swimming systems [Tol 12 1: a central well and bi-even 
symmetry in the x b function, and odd symmetry in the 
y b and 9 b functions. Comparing them with equivalent 
plots for the low Re system [HJ (Fig. [3]B) highlights the 
differences between the environments. Most significantly, 
the x b function's magnitude is 1.5 — 2 times larger for the 
granular swimmer than for the low Re system, in accord 
with the relatively larger F± in the granular medium. 
Interestingly, this relationship is reversed for the 9 b func- 
tions, indicating that turning gaits produce less rotation 
for swimmers in CM relative to those in low Re fluids. 



We note that the zero-sets we identify on the CCFs rep- 
resent the optimal solutions discovered by the parametric 
optimization methods in 
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In summary, CCFs advance our understanding of lo- 
comotion by substituting geometric insight for laborious 
calculation. Since the technique requires only empirical 
force laws, we argue that this method lays the ground- 
work for geometric analysis of biological and robotic loco- 
motion in environments that are not yet (and may never 
be) described by comprehensive equations of motion. We 
propose that the geometric insight gained by CCFs for 
locomotion will be analogous to insight into complex dy- 
namical systems provided by low dimensional maps [27j . 
Future investigations will advance the linear model to 
represent dynamics near geometric singularities and in- 
corporate effort-based distance metrics [IH to describe 
locomotive efficiency. It will also be interesting to com- 
pare optimal locomotor strategies in different media com- 
bining the CCF analysis with effort metrics to identify 
maximally-efficient gaits. Finally, the use of continuous 
curvature modes 13[ will allow application of our frame- 
work to continuous systems like the biological sandfish, 
snakes that move within yielding substrates like loose 
soils, and even nematodes and spermatoza in complex 
biofluids (iij . 
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